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Abstract — Applying the theory of Compressive Sensing in 
practice always takes different kinds of perturbations into con- 
sideration. In this paper, the recovery performance of greedy pur- 
suits with replacement is analyzed when both the measurement 
vector and the sensing matrix are contaminated with additive 
perturbations. Specifically, greedy pursuits with replacement 
include CoSaMP, SP and IHT algorithms, where the support 
estimation is evaluated and replaced in each iteration. Based 
on restricted Isometry property, a unified form of the error 
bounds of these recovery algorithms is derived under general 
perturbations for compressible signals. Derived from the results, 
the recovery performance is stable against both perturbations. 
Also, these bounds are compared with that of oracle recovery — 
least squares solution with the locations of K largest entries in 
magnitude known a priori. The comparison shows that these 
bounds only differ in coefficients from the lower bound of oracle 
recovery with some certain signal and perturbations, as reveals 
that near-oracle recovery performance of greedy pursuits with 
replacement is guaranteed. Numerical simulations are performed 
to verify the conclusions. 

Index Terms — Compressive sensing, sparse recovery, general 
perturbations, performance analysis, restricted Isometry prop- 
erty, greedy pursuits, compressive sampling matching pursuit, 
subspace pursuit, iterative hard thresholding, oracle recovery. 



I. Introduction 

COMPRESSIVE sensing, or compressive sampling (CS) 
|l |-|3|, is a novel signal processing technique proposed 
to effectively sample and compress sparse signals, i.e. signals 
that can be represented by few significant coefficients over 
certain basis. One of the essential issues of CS theory lies in 
recovering A^ x 1 vector x from its linear observations, 



Ax, 



(1) 



where y is an M x 1 vector with M < TV, A is an M x A^ 
sensing matrix, and x is assumed to be iiT-sparse, which means 
only K out of its N entries are nonzero. Unfortunately, directly 
finding the sparsest solution to ([T]i is NP-hard, which is not 
practical for sparse recovery. This leads to one of the major 
aspects of CS theory — designing effective recovery algorithms 
with low computational complexity and fine recovery perfor- 
mance. 
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A. Overview of Sparse Recovery Algorithms 

A family of convex relaxation algorithms for sparse recov- 
ery \Aj had been introduced before the theory of CS was 
established. Based on linear programming (LP) techniques, 
it is shown that ii norm optimization, also known as basis 
pursuit (BP), 



s.t. y = Ax 



(2) 



yields the sparest solution as long as A satisfies the restricted 
isometry property (RIP) with a constant parameter |[5)-||7). 
Recovery algorithms based on convex optimization include 
interior-point methods 18), ||9l and homotopy methods pO) , 
fill. 

In contrast with convex relaxation algorithms, non-convex 
optimization algorithms solve (fill by minimizing £p norm with 
respect to < p < 1, which is not convex. Typical algorithms 
include focal underdetermined system solver (FOCUSS) |12|, 
iteratively reweighted least squares (IRLS) fTSl, smoothed 
£'^ (SLO) |14|, and zero-point attracting projection (ZAP) 
|,15,|. Compared with convex relaxation algorithms, theoretical 
analysis based on RIP shows that fewer measurements are 
required for exact reconstruction by non-convex optimization 
methods |16|. 

A family of iterative greedy algorithms has received much 
attention due to their simple implementation and low computa- 
tional complexity. The basic idea underlying these algorithms 
is to iteratively estimate the support set of the unknown sparse 
signal. In each iteration, one or more indices are added to 
the support estimation by correlating the columns of A with 
the regularized measurement vector. Typical examples include 
orthogonal matching pursuit (OMP) |17|, regularized OMP 
(ROMP) 1 181, and stage-wise OMP (StOMP) |19|. Compared 
with convex relaxation algorithms, matching pursuits need 
more number of measurements, but they tend to be more 
computationally efficient. 

Recently, several greedy pursuits including compressive 
sampling matching pursuit (CoSaMP) |20| and subspace pur- 
suit (SP) pT) have been proposed by incorporating the idea 
of backtracking. In each iteration, SP algorithm refines the K 
columns of matrix A that span the subspace where measure- 
ment vector y lies. Specifically, SP adds K more indices to 
the K candidates, and remaining the most reliable K ones. 
Similarly, CoSaMP adds 2K more indices in each iteration, 
while computing the regularized measurement vector in a 
different way. By evaluating the reliability of all candidates 



in each iteration, these algorithms can provide comparable 
performance to convex relaxation algorithms, and exhibit low 
computational complexity as matching pursuit algorithms. 

Another kind of greedy pursuits, including iterative hard 
thresholding (IHT) fS^l and its normalized variation f23l, is 
proposed with the advantages of low computational complex- 
ity and theoretical performance guarantee. In each iteration, 
the entries of the updated solution except for the most K 
reliable ones are set to zero. Together with CoSaMP and SP, 
these algorithms can be considered as greedy pursuits with 
replacement involved in each iteration. 

B. Overview of Perturbation Analysis 

To apply the theory of CS in practice, the effect of noise 
must be taken into consideration. The most common analysis 
is the additive noise to the measurements, i.e. 



y = y 



(3) 



where e is termed measurement perturbation. Most existing 
algorithms have been proved to be stable from perturbed 
measurements, including theoretical analysis of BP (|6), ROMP 
1|24), SP 1|2T1, CoSaMP [W], and IHT |22), ||23|. It is shown 
that the error bounds of the recovered solutions are propor- 
tional to the energy of measurement perturbation e. A certain 
distribution of e can be introduced to achieve better results, 
such as Gaussian |,25J-[28J or others [29 J , [30J . 

Until recently, only a few researches involve the perturbation 
to A, which is also termed system perturbation. Representative 
researches include analysis of BP [311, CoSaMP ||^, SP | 
ip norm optimization |[34) with < p < 1, and OMP | 



In their studies, ([3]l is extended by introducing a perturbed 
sensing matrix, i.e. A = A + E. It is of great significance 
to analyze the stability of recovery algorithms against both 
perturbations when the theory of CS is applied in practice. 
Other related works include mismatch of sparsity basis [ ,36| 
and sparsity-cognizant total least-squares p7) . 

Two practical scenarios are usually considered in CS appli- 
cations. First, when A represents a system model, E denotes 
the precision error when the system is physically implemented. 
Thus the whole sensing process is 



y = (A + E)x- 



(4) 



and only the nominal sensing matrix and contaminated mea- 
surement vector are available for recovery, i.e. x = R(A,y), 
where R(-) denotes a certain recovery algorithm. 

In countless other problems, E is involved due to mis- 
modeling of the system A. Thus A = A + E and the sensing 
process is 



y = Ax + e. 



(5) 



Both the sensing matrix and measurement vector are contam- 
inated, and the recovered solution is x = R(A,y). 

Available researches ||3T|-||33l, |J35| are all based on the 
latter scenario, thus dSll is well considered and fully analyzed 
in this paper. The first scenario (pi is briefly discussed after 
that as a remark. 



C. Our Work 

This paper mainly considers the recovery performance of 
greedy pursuits with replacement against general perturba- 
tions. Specifically, when both measurement perturbation and 
system perturbation exist, the error bounds of the solutions of 
CoSaMP, SP and IHT are derived and discussed in detail for 
not strictly sparse signals, i.e. compressible signals. It is shown 
that the bounds of these greedy pursuits are proportional to the 
relative perturbations of both sensing matrix and measurement 
vector, which indicates that the recovery performance is stable 
against general perturbations. These error bounds are also 
compared with that of oracle recovery, and it is shown that 
the results in this paper are optimal up to the coefficients. 
Numerical simulations are performed to verify the conclusions 
in this paper Previous and related works including p8) , p2) , 
| |33| are introduced in Section VI. 

The remainder of this paper is organized as follows. Section 
II gives a brief review of RIP and three greedy pursuits 
with replacement: CoSaMP, SP, and IHT. Section III presents 
the main theorems when greedy pursuits with replacement 
meet general perturbations, and the results are shown to be 
near-oracle by comparing with oracle recovery in Section IV. 
Numerical simulations are performed in Section V to verify the 
contributions in this paper Several related works are discussed 
in Section VI, and the paper is concluded in Section VII. 

II. Preliminary 

In this section, the definition of RIP and descriptions of 
several greedy pursuits, including CoSaMP, SP, and IHT, 
as well as their recovery performance against measurement 
perturbation, are introduced. 

A. Restricted Isometry Property (RIP) 

The restricted isometry property (RIP) for any matrix A 
describes the degree of orthogonality among its different 
columns. 

Definition 1 For each integer K = 1,2,..., define the 
restricted isometry constant (RIC) 5k of a matrix A as the 
smallest non-negative number such that 



(1-5k)||x||^<||Ax||^<(1 + <5k)||x||^ 
holds for all K-sparse vectors x. 



(6) 



According to Definition I, 5k < 1 implies that every K 
columns of A are linearly independent, and 52k ^ 1 implies 
that A almost maintains the £2 distance between any pair of 
iiT-sparse signals. It is also easy to check that if A satisfies 
the RIP with 5ki and Sk2, and Ki < K2, then Sk^ < 5k2- 

Calculating the exact value of RIC given a fixed matrix 
A is intractable because it involves all submatrices formed 
by taking K columns of A. Fortunately, random sensing 
matrices have small RICs for overwhelming probability. For 
example, if the entries of A e M^^^^ are independently and 
identically distributed Gaussian random variables of zero mean 
and variance l/M, and M satisfies 

,, CKlog(N/K) 



then 6k < £ with at least probability 1 — c ^^^, where C and 
c are two constants 15|, | |20l . Recent results about the bounds 
of RIC can be found in ||38|, (|39). 

Though RIP is usually used to show the impact of sensing 
matrix on sparse signals, the following lemma permits us to 
generalize the results for sparse signals to general signals. 



TABLE I 

The Specification of the Constants 



CoSaMP 



SP 



1.26 



IHT 



Lemma 1 (Proposition 3.5 in [20^) Suppose the sensing 
matrix A satisfies RIP of sparsity level K with 5k, then for 
every signal x G M^, it holds that 



iiAx|i2< v/m^fiix 



1 



K 



(7) 



Lemma [T] is used in the following analysis for compressible 
signals. 

B. Greedy Pursuits with Replacement 

Greedy pursuits with replacement include CoSaMP, SP, and 
IHT algorithms, and they are briefly introduced as follows. 

The CoSaMP algorithm iteratively refines the support of K- 
sparse vector x. In each iteration of CoSaMP, 2K more indices 
are selected by correlating A with residue vector r, then the 
best K candidates out of at most "iK ones are remained and 
the residue vector is updated. The details of CoSaMP can be 
found in f20l. 

The SP algorithm is firstly proposed in [40] and further 
developed in pTj . It iteratively refines the K columns that 
span the subspace where the measurement vector y belongs. 
Unlike CoSaMP, SP adds K more indices in each iteration, 
and remains the best K candidates. In addition, the regularized 
measurement vector of SP is orthogonal to the subspace 
spanned by the columns of A indexed by the K candidates, 
making the new K indices added in each iteration totally 
different from previously identified K ones. Please refer to 
[ |2T| for more details. 

The IHT algorithm is firstly proposed in pT) , and later 
developed and analyzed in fTT\. To guarantee the convergence 
of the method, a normalized variation NIHT is proposed in 
[ [23] which retains theoretical performance guarantee similar 
to IHT. Without loss of generality, only IHT algorithm is 
discussed in the following. 

The following Theorem [T] reveals the error bound of the 
solution of greedy pursuits with replacement when only mea- 
surement perturbation exists. 

Theorem 1 Given a perturbed measurement vector y = 
Ax + e where x is K-sparse, the estimated solution x^'l in 
the l-th iteration of CoSaMP and IHT algorithms satisfies 



|x-xW| 



<C||x-x['-il| 



Cillell 



(8) 



and the estimated support set T' in the l-th iteration of SP 
algorithm satisfies 



|XT-T'll2 < C'l|xT-T'-i|l2 + C*! 11*^11 



(9) 



where T denotes the support of x, i.e. the set of locations 
of nonzero entries of x, T — T' denotes the set comprised of 



C 



Ci 



D 



0.171 
4(54jf 



(1 - SiKV 

6 + 2S4K 
(1 - S^kY 



0.206 
2<53K + 2(53^ 

(1 - S3k)'' 

4(1 + <53k) 
(1 - SiK^ 



1_6,53K+6«|~ 



■^^Ik+^Ik 



0.353 

VSSsK 

2^1 + Sk 

2 + S3K 

i-VaSsK 



a; e T and x ^ T', and xt denotes the subvector composed 
of entries of x indexed by set T. 

Furthermore, if the sensing matrix A satisfies dbK < c, then 
C < 1, and it can be derived that 



|x-xW| 



<aC" llxll 



DM 



(10) 



holds for greedy pursuits with replacement. The specific values 
of the constants a, b, c, C < 1, Ci and D are illustrated in 
TABLED 

Proof: For CoSaMP algorithm, inequality ([8]l can be 
obtained by following the steps of proof of Theorem 4.1 in 
|20|, except for preserving RICs during the derivation. By 
recursion, inequality (10 1 can be proved from ([8]l. For SP 
algorithm, the result cannot be directly found in 121]. However, 
it can be proved based on several facts about SP algorithm 
in pT) , and the proof is provided in Appendix A. For IHT 
algorithm, inequality (|8]lis derived in the proof of Theorem 5 
in |22J, and inequality ( 10 1 can be proved by recursion. ■ 



C. Discussion 

Now we have introduced three greedy pursuits: CoSaMP, 
SP, and IHT, where sparsity level K needs to be known a priori 



for replacement. According to ( 10 1, it can be concluded that 



when only measurement perturbation exists, the error bounds 
of them are exponential decay function of iteration number, 
and the limit of these error bounds is proportional to the energy 



of perturbation. Two remarks based on (lOi are derived as 
follows, and the main conclusions in this paper are based on 
( fTO] ) as weU. 

Remark 1 When no perturbation exists, then after finite 
iterations, the recovery solutions of CoSaMP and SP are 
guaranteed to be identical to the original signal x. The result 
can be verified through the following statement. 

Suppose .Tmin IS the smallest magnitude of the nonzero 
entries of x. Then after 



l^ 



logc 



1 



iterations, the recovery solution xW obeys 



< x„ 



(11) 



If the support is perfectly recovered, then for CoSaMP and 
SP, the solution is akeady identical to the original signal. 
Otherwise, at least one nonzero entry is not detected, thus 
the recovery error is no less than Xmin, which contradicts 
^n\ . Notice that the solution of IHT does not have the above 
property, since exact support recovery does not imply exact 
signal recovery. 



Remark 2 According to (10 1, after 



/ = 



logc 



.a||x||2. 
iterations, the error bound of the recovery solution satisfies 

||x-xW||2<(i? + l)||e||2, (12) 

which means that error bound is proportional to the energy of 
perturbation, and the recovery performance of greedy pursuits 
with replacement is stable against measurement perturbation. 
It needs to be emphasized that the above analysis is per- 
formed as worst-case, and the demands of RICs of sensing 
matrix are sufficient conditions. Better recovery performance 
of the three algorithms is normally achieved in practice, thus 
they are widely used in various applications. 

III. Recovery Performance of Greedy Pursuits 
WITH Replacement 

In this section, the error bounds of greedy pursuits with 
replacement are derived for general signals when both mea- 
surement perturbation and system perturbation exist. 

A. Notations and Assumptions 

Several notations and assumptions are stated first and they 
will be used in the following analysis. The original signal x 
is assumed to be compressible signal in the following context, 
which means that it can be well approximated by a sparse 
signal. Xx denotes the subvector composed of entries of x 
indexed by set T. Ax denotes the submatrix composed of 
the columns of A indexed by set T. A^^ denotes the pseudo- 
inverse of a tall, full-rank matrix A by A^ = (A*A)^^A*. 
Vector xk is assumed to be the best K-term approximation 
to X, and define the approximation error x^ = x — xa'. The 
approximation error can be quantified as 



rK 



^Kll2 



SA' 



•-KWl 



(13) 



When the original signal x is X-sparse, the ratios tk = sk — 
0. If X is compressible, then for reasonable K, the ratios tk 
and Sk are expected to be far less than 1. 

The symbol |jA||2 denotes the spectral norm of a matrix 
A, and the symbol ||A|J2 denotes the largest spectral norm 
taken over all iiT-column submatrices of A. The measurement 
perturbation e and system perturbation E can be quantified 
with the following relative bounds: 



<ev 



lEl 



where 



|y||2, 



< sa, 



\{K) 



lEI 



m 



<£ 



(K) 



(14) 



7^ 0. Ey quantifies the relative 
measurement perturbation, while ea and e)^ ' quantifies the 



relative system perturbation. In practical scenarios, the exact 
forms of e and E are not known in advance, thus the relative 
upper bounds are applied instead. In this paper, Cy, e^^, e)^ 
are assumed less than 1. 

The following lemma quantifies the RIC of the perturbed 
sensing matrix A = A + E. 



Lemma 2 (Theorem 1 in [31^) Assume the sensing matrix A 
satisfies RIP of sparsity level K with Sk, and the relative 
perturbation e)^ is associated with matrix E. Define the 
constant 



OK,, 



^(1 + <5a) 1+4 



J-Fc; 



1, 



(15) 



then the RIC Sj^ of matrix A = A + E satisfies 5k < (^A-^max- 

As can be seen from Lemma l2] the upper bound of the 
RIC of perturbed matrix A is slightly bigger than that of A 
if the relative perturbation ej^ is small. Notice that ^A.max 
represents a worst-case of 5k, thus better RIC of A is normally 
achieved in practice. 

B. Error Bounds of Greedy Pursuits with Replacement 

The relative error bounds of greedy pursuits with replace- 
ment under general perturbations for compressible signals are 
given in this subsection. The following Theorem [2] summarizes 
the main result. 

Theorem 2 Suppose y = Ax, where x is an arbitrary vector, 
and the available information used to recover x k y = y + e 
and A = A + E, where e and E denote measurement pertur- 
bation and system perturbation, respectively. If the available 
perturbed sensing matrix A satisfies RIP with 



5bK < c, 



(16) 



and the non-perturbed matrix A satisfies RIP of order K with 
5k, then in the l-th iteration, the relative error of the solution 
xl'l of greedy pursuits with replacement satisfies 

llx-xW|| 



TK 



+ {l + Sy){rK + SK) 



K 



AK) 



where the specific values of the constants a, b, c, 
are illustrated in TABLE U\ and TABLE Wl\ respectively. 
Furthermore, after at most 



(17) 
C <l, D 



IK) 



logc 



SA 



(18) 



iterations, the greedy pursuits estimate li. with accuracy 



rWl 



< 



b^ 



TK 



-{K 



+ 4 + (1 + £y) (^'A' + sk) 



(19) 



For exact sparse signals, a better and more intuitive result 
can be derived from Theorem l2] by setting vk and sk to zero. 
The result is stated as Theorem |3] 



TABLE II 

The Specification of the Constants 



CoSaMP 



SP 



IHT 



C 



4(54 



2<53A' 



(1 - hK? 
£) 6+2<54jf 



2^|k 



(1 - hKY 
l-6i3K+64?R^-253 - 



V853K 
2+^3Jf 



Theorem 3 Suppose y = Ax, where x is K-sparse, and the 
available information used to recover yi is y = y + e and 
A = A + E. /f the available perturbed sensing matrix A 
satisfies RIP with 



ObK 



<c, 



and the non-perturbed matrix A satisfies RIP of order K with 
5k, then in the l-th iteration, the relative error of the solution 
xl'l of greedy pursuits with replacement satisfies 



|x-xW| 



< a& + Dy/l + SK (sy + e^ 



iK)\ 



(21) 



where the specific values of the constants a, b, c, C < 1, and 
D are illustrated in TABLE I//I 
Furthermore, after at most 



l^ 



JK)' 



logc 



(22) 



iterations, the greedy pursuits estimate x with accuracy 



rW\ 



< 



{dVi 



3 A' 



(e. + 4"0 



(23) 



C. Proof of Theorem [2] 

Proof Recalling the sensing process dSll with general 
perturbations, it is equivalent to 



y = Ax + e = f A - El (xa- + x.%^) + e 

= AxA' + (e-ExA + Ax^). (24) 

Define e = e — Exa + Ax^ as the error term. According to 
(10 1, under certain condition (16i about the RlCs of A, the 



solution xW in the l-th iteration satisfies 



XA 



rWl 



<aC' ||xA||2 + £'||e| 



2 ' 



(25) 



where C < 1 and D are constants specified in TABLE |II 
From the triangle inequality, it can be derived that 



rWI 



< IIx-xaI 



xa -X 



[i]\ 



Substituting ( [26] l into ( [25] l, and from the fact that 

IIxaIIs < l|x||2, 

the error bound in the l-th iteration satisfies 



rWl 



<aC'||x|| 



^Kll2 + ^l|e|l2- 



(26) 



(27) 



(28) 



To estimate ||e||2, the triangle inequality. Lemma [T] and 
definition (13i implies 



||e||2<||e|l2 + |lExA||2 + ||Ax^.||2 



<lle|| 



lEllf^llxAl 



W 



||x||2(rA + SA). 



(29) 



It can be derived that 



|e||2 <£y ||y|l2 = £y l|Ax||2 < Ey (||Axa||2 + l|Ax^||2) 



<£y||A||f)||x||2(l + rA + SA), 



and 



2 



'A 



\iK) 

2 



(30) 



(31) 



and (28 i, the inequality ( 17 1 can be obtained 



(20) Substituting (|30]) and (mi into ([29]l, and together with p7| ) 



For the second part of the theorem, notice that when ( 1 8 1 
holds, it is obvious that 



aC' 



<£v 



+ e^^ + (l + £y)(rA + SA) 



(32) 



and ( [T9| follows immediately, which completes the proof of 
Theorem |2l ■ 



D. Discussion 

In this subsection, several remarks are drawn to have a deep 
insight into the results of this paper. 

Remark 3 As can be seen from ([T9|, the relative error of 
greedy pursuits with replacement is proportional to both kinds 
of perturbations and approximation error of x to iiT-sparse sig- 
nal. If the original signal is j^-sparse and the sensing process 
is non-perturbed, according to (21 1 with £y = e)^ = 0, the 
relative error is bounded by an exponential decay function of 
iteration number, which indicates that the recovery solution 
can approach the original signal at any given precision. For 
exact iiT-sparse signal x sensed under general perturbations, 
( [23] ) demonstrates that the recovery accuracy is limited by 
the size of both perturbations, and the recovery performance 
is stable against general perturbations. If the sensing process 



is non-perturbed for general signal x, according to (19i with 
Ey = £)^ = Q, the recovery accuracy is determined by how 
well it can be approximated with A' -sparse signal. 

Remark 4 Under the assumption ( 16 1, the base C of the 
exponential function in ( [T7| ) is less than one, which guarantees 
the convergence of the recovery error. As ^^a" — ^ c, the base 
C approaches one and the constant D approaches infinity. 
Also, due to the monotonicity of RIC, for less sparse signal, 
i.e. larger K, the constant D gets larger, which means the 
proportionality constant is larger. 

For these three greedy pursuits with replacement, K is quite 
a significant parameter and needs to be known a priori. An 
improper selection of K may result in great loss of their 
recovery performance. If K is selected too large, the condition 
( 16 1 can hardly be satisfied, and the parameter D may be quite 
large. On the other hand, too small K indicates large tk and 
sk, which also leads to poor performance. 



As a typical compressible signal, strong-decaying signal x 
is a vector whose ordered coefficients satisfy 



>p|x|, 



KK N - 1, 



(33) 



where |xLn is the ^th largest magnitude among the elements 
of X, and p > 1 is a constant controls the speed of the decay: 
the larger p is, the faster x decays. Simple calculation implies 
that 



^ -K ^ Vp+^ -k 

For compressible signal satisfies a power law, i.e. 



1(0 



R-l- 



1<1<N, 



(34) 



(35) 



where i? is a positive constant denoting the radius of weak- 
£i/p ball, and p > 1 controls the speed of the decay. It can 
also be calculated that 



rK^K'/^-P, SK- "^^^ / j^i/^'P. (36) 

p-1 

In fact, the parameter K should be adjusted according to 
£y, Ejs^ and the decay parameter p. In the scenario with 
relatively large perturbations, according to ( [T9] l, K can be 
selected a small value such that S^x and D are small, and 
rji and Sk are comparable with Sy and e^ \ In the scenario 
with relatively small general perturbations, the recovery error 
is mainly determined by vk and sk, and K is a trade-off 
parameter between these ratios and D. 

Remark 5 Considering the sensing process (|4|) with A = 
A + E, 



y=(A + E)x- 



e = AxA' + ( e + ExA- + Ax^ 1 . (37) 



Define e = e + Exa' + Ax^ as the error term. Following 
the steps of proof of Theorem l2] it can be derived that if the 
sensing matrix A satisfies RIP with SbK < c, then 

rWll <aCM|x|L 



X — x^ 

M II z 

Since 

Il§ll2<l|e|l2 + I|E|,; 
and the fact that 

||-f*-||2 — 11-^112 



'■A'112 



D e 



(38) 



,(K) 



X7< 



|A||f^||x||2(rA + SA) 



lEI 



W 



< 1 



JK) 



\{K) 



(39) 



it can be derived that the relative error of the solution in the 
l-i\\ iteration obeys 



Ix-xWl 






IK) 



(40) 

Notice that in this scenario the error bounds of greedy 
pursuits with replacement are also proportional to the relative 
perturbations of both sensing matrix and measurement vector, 
which indicates that the recovery performance is stable against 
general perturbations. Also, it should be stressed that for both 
scenarios, the demands of RIP are all for the available sensing 
matrix in the process of recovery. 



IV. Comparison with Oracle Recovery 

In this section, the upper and lower bound of oracle recovery 
is derived, and is compared with those of greedy pursuits with 
replacement. Then it is shown that the results in this paper are 
optimal up to the coefficients. 



A. Error Bound of Oracle Recovery 

Consider the oracle recovery where the locations of the K 
largest entries in magnitude of a vector x are known a priori. 
Assume T is the set of the K locations, i.e. T = supp(xa). 

Recall the sensing process (|5]l where both perturbations ex- 
ist. Through least squares (LS) method, the estimated solution 
X of oracle recovery is obtained by 



;t - 



XT = A^y, XT- = 0. 
It is easy to check that 



XT 



A^y 



:t 



AxR- 



XT + ALe. 



Thus the estimation error of the LS solution obeys 



l|x-x|L = ||XT 



- I|2 
XTII2 



I Xt -l|2 



f IIXTo 
^k\\2 ■ 



- I|2 



According to (43 1, the estimation error obeys 



< 



< 



|A^|l2l|e|l2 



*-_R-|l2 



'-K\\2 ■ 



Substituting (|29j into (|44]) and from the fact that 

1 



lA^II < 



Vi-s 



= D, 



(41) 



(42) 



(43) 



(44) 



(45) 



K 



the relative error of LS solution is derived following the steps 
in the proof of Theorem [2] 



<{d^/1 + 8k + 1 



+ e^^ + (1 + Ey) {rK + sk) 



(46) 



Comparing ( |46] l with ( [T9| , it can be derived that after 
finite iterations, the error bounds of greedy pursuits with 
replacement and the error bound of oracle recovery only differ 
in coefficients. The error bounds of LS solution and greedy 
pursuits with replacement are all on the same order of the 
noise level and approximation error 



0(£y + £f' + (l + £y)(rA + .SA)) 

under general perturbations. 

B. Lower Bound Analysis of Oracle Recovery 



(47) 



According to (|28]l and (|29]l, after 



l^ 



log 



6 2 



IE 



W 



XK 2 



■c 



a||A||f)||x||2 



iterations, the recovery error of greedy pursuits with replace- 
ment satisfies 



rWl 



< D 



1 






lEllf^llxxl 



(48) 



To show the optimaHty of the results in this paper, we prove 
that for some e, E and x, the lower bound of estimation error 
of oracle recovery is also dependent on the following three 
terms 



lEI 



(A') 



|xk| 



||A||f)||x||2(r,^ + ,SK), (49) 



which denote measurement perturbation, system perturbation 
and approximation error, respectively. 

According to (|43]l, the estimation error obeys 



(50) 



■x||2 > \\A\,e\\^ > 



where the second inequality is due to the fact that 

1 



A. 



Atl> 



where Amin denotes the smallest singular value. 

In particular, the system perturbation E can be chosen such 
that 



||ExK|l2 = |lE||f)||xK|l2, 
Meanwhile, it can be assumed that 

(-Exk)* (Ax^) > 0, 



(51) 



(52) 



otherwise — E is applied and pT) still holds. Furthermore, the 
measurement perturbation e can be chosen to be orthogonal 
to Ex/< and AxJ^. Due to the above assumptions, it can be 
obtained that 



> 
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(53) 



For the case that xj^ is iiT-sparse, according to Cauchy- 
Schwartz inequality 



> 



II^A'll2 

it is derived that 
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(55) 



Substituting (51 1 and (55i into (53 i, and according to (50i 



the estimation error of oracle recovery is lower bounded by 

I2>^lfl|e|l: 
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(56) 
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Fig. 1. Average relative error (1000 trials) of greedy pursuits with replace- 
ment and oracle LS recovery versus sparsity level K when ea. = 0.05 and 
Ev = 0.05. The RIP condition seems to be satisfied in this scenario. 



where 



Di^ 



1 



%/r^^ 



(57) 



Comparing (J56]) with (J48]l, it can be concluded that there 
exist X, e, and E for which the results in this paper are 
essentially tight up to the coefficients. In general, no recovery 
algorithms can do better than the oracle least squares method. 
Thus, the greedy pursuits with replacement: CoSaMP, SP, and 
IHT, can provide near-oracle recovery performance against 
general perturbations. 

V. Numerical Simulations 

Three numerical simulations are performed in MATLAB 
and demonstrated as follows. The first simulation compares 
the relative error of greedy pursuits with replacement and 
oracle LS method for different sparsity levels. The second 
simulation compares the recovery error for different relative 
perturbations. The final simulation considers compressible 
signal and presents the influence of parameter K on recovery 
error. 

In each trial a new matrix A of size 512 x 2048 is randomly 
generated with independent Gaussian distributed entries with 
cr^ = 1/512 so that the expectation of the I2 norm of each 
column is normalized. The nonzero locations of a sparse 
vector X are randomly selected among all possible choices, and 
the nonzero entries are generated independently from normal 
distribution. Then, for each pair of relative perturbations, A 
and y are generated according to (fSl) and (fT4|. Notice that 

— (K) 

only ea is used in the simulations since calculating ej^ ' is 
NP-hard. The replacement is reasonable because Ea ~ £a 
holds for all K with high probability when both A and E are 
random Gaussian matrices | |31J . 

In the first experiment, the relative perturbations are set to 
£a = 0.05 and Sy = 0.05. The sparsity level K varies from 
1 to 25, and the simulation is repeated for 1000 times. Fig. [T] 
presents the average relative error for the range where the 
demand for RIP seems to hold. As can be seen from this figure. 
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Fig. 2. Average relative error (1000 trials) of greedy pursuits with replacement and oracle LS recovery versus different £a and Sy. The results are represented 
for the same K in the same row, and for the same recovery algorithm in the same column, with labels in the left and on the top, respectively. 



when K < 20, the recovery error increases almost linearly as 
K increases. When the sparsity level K continues to increase, 
the recovery performance of IHT cannot be bounded any more, 
which indicates that the RIP condition is violated. 

In the second experiment, for fixed sparsity level K — 5, 
10 and 15, the relative perturbations ea and Sy vary from 
to 0.1 with step size 0.005, respectively. The simulation 
is conducted 1000 trials to obtain the average relative error 
The results are demonstrated in Fig. [2J where the results are 
represented for the same K in the same row, and for the 
same recovery algorithm in the same column, with labels in 
the left and on the top, respectively. As can be seen from 
them, for the same K, the surface of oracle LS method lies 
below the others, while the surface of CoSaMP is above on 
the top. The surfaces of SP and IHT are almost the same in 
the middle. It verifies that the relative errors of CoSaMP, SP, 
and IHT are proportional to the measurement perturbation as 
well as system perturbation, and these greedy pursuits with 
replacement can provide near-oracle recovery performance. It 
needs to be pointed out that for the same recovery algorithm, 
the gradients of these surfaces increase with the increase of 
K, which confirms that the coefficients of error bounds are 



related to K as discussed in Section III and IV. To display the 
results more intuitively, the average relative error is plotted 
versus ea for different K when Sy = 0, as in Fig. l3] As can 
be seen, the recovery error scales almost linearly as Ea, and 
the slopes of these curves increase as K for the same recovery 
algorithm. 

In the final experiment, we consider compressible signal x 
satisfying a power law 



1(0 



= 1- 



l<l< N, 



(58) 



where p — 2 in this experiment. The relative perturbations are 
set to Ea — 0.01 and £y = 0.01. The parameter K in the 
recovery algorithms varies from 5 to 100, and the results are 
plotted in Fig. [4] Since IHT fails when K > 20, the recovery 
error of it is not shown for K > 20 in the figure. As can be 
seen from Fig. p] the optimal K is around 22 for CoSaMP 
and SP, 19 for IHT, and around 40 for oracle recovery. Too 
small K or large K will affect the recovery performance for 
compressible signals, as is discussed in Remark 4. According 
to (36 1, for K = 22 and p = 2, the approximation error r^- 



is about 0.010, and sk is about 0.017, which are comparable 
with Ea and Ey. 
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Fig. 3. Average relative error (1000 trials) of greedy pursuits with replace- 
ment and oracle LS recovery versus £a for different sparsity level K when 
£v =0. 
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Fig. 4. Relative error of greedy pursuits with replacement and oracle LS 
recovery for compressible signal versus different parameter K when ea = 
0.01 and Cy = 0.01. IHT algorithm becomes invalid when K > 20, thus it 
is not plotted for the wider range. 



VI. Related Work 

A similar concept of this paper is presented in 1*281, where 
the authors also establish a near-oracle performance guarantee 
for CoSaMP, SP, and IHT algorithms. Based on RiP the 
analysis considers the recovery of a X-sparse signal with 
the assumption that the measurement vector is corrupted by 
additive random white Gaussian noise. The main result of | |28| 
is stated as follows, where some notations are replaced for the 
sake of consistency. 

Assume that the white Gaussian noise vector is with co- 
variance matrix cr^I, and that the columns of sensing matrix 
A are normalized. Under certain conditions of RICs and with 
probability exceeding 1 — [^Jtt[1 + a) log A^- N°-)^^, it holds 
that 



x-x||^ < C- {2{l + a)\QgN) ■ Ku\ 



(59) 



where x is the recovered signal, and C is a constant related 



to recovery algorithm, sensing matrix and sparsity level. The 
result is similar to those for the dantzig-selector and the basis 



pursuit, but with different constants. The logA^ factor in (59 1 
is proven to be unavoidable in |42|, therefore this bound is 
optimal up to a constant factor. The result is also extended to 



the nearly-sparse case in |28|. 

Other relevant researches are presented in p2| , p3) , where 
the error bounds of CoSaMP and SP are also derived under 
general perturbations. The analysis is performed following the 
steps of 1 3 1 1 for basis pursuit. Define 
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X-X^ 2 
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X-X/f 1 
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||XK||2 ■ ^/K\\^K\W %/r^^' 

The main result in 1*321 states that under the conditions that 
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'3AK 
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and that 



ax + (3k < 
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the recovered solution of CoSaMP satisfies 

1 ,, 
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IIEI 



(60) 



(61) 



||2<C- ||X-XA-||2 
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\2 
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Further define 



where 



SK)(K) 




£y lly||2, 



(62) 



(63) 



(64) 



The main result in |l33l states that under the conditions that 

1.083 



33a: 



< 



(i+^r'F 



1 



and that 



aK + Pk < 



AK)-' 



(65) 



(66) 



the recovered solution of SP satisfies 

||x-x||2 <(Ca + 1)||x-xa||2 + Ca 

+ Gi^ep^^K,y 



X-XK 1 



K 



(67) 



Our analysis, on the other hand, derives a unified form of 
the error bounds of CoSaMP, SP, and IHT under general per- 
turbations. Specified in Theorem |2] for compressible signals, 
when measurement vector and sensing matrix are perturbed, 
under the condition of ( [T6] l, the error bound of recovered 
solution in each iteration is derived as (17i. It is also proved 



that after finite iterations, the error bound is proportional to 
both perturbations and approximation error, as in ([T9]l. The 
result is compared with oracle recovery, where the locations 
of K largest entries in magnitude are known a priori. It is also 
proved that for some certain x, e, and E, the lower bound 
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of oracle recovery only differs in coefficients from the error 
bounds of greedy pursuits with replacement. Therefore, the 
near-oracle recovery performance of them is guaranteed in this 
work. 

The main difference between (28) and this work is that we 
consider a more general completely perturbed scenario, and the 
optimality of the recovery performance is also in this sense. 
Compared with p2| , (331, the demand of RIP in our work 
is for the perturbed sensing matrix A other than A, which is 
due to the fact that only A is available for recovery. Also, 
the demand of RICs is with a constant parameter here. In 



addition, the condition such as (61 1 or (66 1 is not required in 
our assumption. Our result is compared with oracle recovery 
and shown to be optimal up to the coefficients, and it is verified 
by plentiful numerical simulations. 



and therefore 



Lemma HI directly implies inequality (J9]l. It's easy to check 
that if d^K < 0.206, then C < 1, and 



jXrp_rpl||2 < C IIXX — T"! 



2 + Y^I|e|l2 



< C" ||x||2 + ^-^ ||e,i2 



(70) 



Applying Lemma[3]to (|70|l and with the fact that S^k < 0.206, 
inequality ( lOi can be derived. 
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VII. Conclusion 

In this paper, the error bounds of three greedy pursuits 
with replacement under general perturbations are derived. It is 
shown that these bounds are proportional to both perturbations, 
which leads to the stability of these algorithms against general 
perturbations. Furthermore, these error bounds are compared 
with that of oracle LS solution with the locations of K 
largest entries in magnitude known a priori. Since the error 
bounds of CoSaMP, SP, and IHT algorithms are on the same 
order of the lower bound of oracle solution for some certain 
signal and perturbations, it can be concluded that these greedy 
pursuits with replacement can provide near-oracle recovery 
performance against general perturbations. Numerical simula- 
tions verify that the recovery errors of these greedy pursuits 
are proportional to both perturbations, and they exhibit near- 
oracle recovery performance. Discussions and simulations also 
reveal how moderate parameter K achieves better recovery 
performance for compressible signals. 

Appendix A 
Proof of Theorem[T]for SP Algorithm 

Several facts about the SP algorithm are firstly introduced 
and proved in pT[, thus they are stated here without proofs. 



Lemma 3 (Lemma 3 in [21]) Let x G M^ be a K-sparse 
vector, and let y — Ax+e be a perturbed measurement vector 
for which A g jjAfxw ^Qf/^yjgj. j/ig j^jp with parameter 5k- 
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